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We trace the evolution of a spherically symmetric density perturbation 
in the MOdified Newtonian Dynamics (MOND) model. The background 
cosmological model is a A-dominated, low-fib Friedmann model with no 
Cold Dark Matter. We include thermal processes and non-equilibrium 
chemical evolution of the collapsing gas. We find that the first density per- 
turbations which collapse to form luminous objects have mass ~ 10 5 M Q . 
The time of the final collapse of these objects depends mainly on the value 
of the MOND acceleration ao and also on the baryon density £lf For the 
"standard" value ao = 1.2 x 10~ 8 cm/s 2 the collapse starts at redshift 
z ~ 160 for tt b = 0.05 and z ~ 110 for fl b = 0.02. 

PACS numbers: 95.30.Lz, 95.30.Sf, 98.35.Mp, 98.80.Bp 

1. Introduction 

Recent developments in cosmological observations have led to so-called 
cosmological concordance model with about 0.03, Q m (dark+baryonic) 
about 0.3 and f^A about 0.7. However, as the ACDM models are dominated 
by hypothetical vacuum energy and non-baryonic dark matter contributions, 
some scientists look for different solutions. Perhaps the most intereresting 
alternative model is the MOdified Newtonian Dynamics model (MOND) 
proposed by M.Milgrom[l]. It assumes that there is no non-baryonic dark 
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matter (or it is negligible) and the lack of matter is only apparent due to 
modification of dynamics or gravity for small accelerations (a <C ao where 
ao is some constant). This model seems to work very well for spiral galaxies 
and many other types of objects[2] but, however, it has some unresolved 
problems (e.g. lack of covariance) . 

Our aim there is to study what would be the implications of the MOND 
model for the formation of the very first objects in the Universe. 



2. MOND vs the standard theory of linear perturbations 

To apply the MOND model to structure formation calculations one en- 
counters a number of difficulties. First of all, MOND is not a theory, it 
is rather a phenomenological model. In its present form MOND is incon- 
sistent with the General Relativity. Up to now there were no successful 
attempts to find a generally covariant theory that could be a generalisation 
of the General Relativity and would give a MOND-like predictions in the 
low-gravity limit [2]. 

MOND is a model that modifies either dynamics or gravity (in this paper 
we assume this second possibility). It introduces a new fundamental scale, 
usually called ao- Gravitational fields much stronger than ao are identical 
to their Newtonian limit gjy and very weak fields are y/aogN- According to 
Sanders and Verheijen [3] the value of the fundamental acceleration scale is 
ao = 1.2 x 10~ 8 cm/s 2 . More precisely, the strength of the gravitational field 
may be written as 

— ) 9 = 9n, (1) 
\ao/ 

where is some function that interpolates between these two extreme 
cases. This function is not specified in the model. We have decided to apply 
the function used by Sanders and Verheijen [3]: 

Kx) = ~^=^ (2) 



VI + x 2 



and, finally, 



9 = 9N\ - (3) 

where x = /ao. 

If we consider the gravitational field of a point mass M, for distances 
R > y/GM/ao the gravitational field would be in the MOND regime, where 
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its strength would decrease as 1/R instead of l/R 2 . It means that there is 
no escape velocity and all systems are gravitationally bound. 

The conseqences of the MOND for cosmology are not studied in de- 
tails yet. R.H.Sanders [4] suggested that because in the early Universe the 
MOND radius is much lower than the radius of the horizon the evolution of 
the scale factor is described by the standard Friedmann equations. Here we 
follow this assumption and study the formation of the first objects in the 
Universe with modified dynamics. 



Let us consider a homogenous ball of density g and some radius R, 
expanding uniformly in all directions with speed proportional to the distance 
from the center. For a sphere of radius r the deceleration in the Newtonian 
limit is 



and, of course, expansion is scale-invariant because deceleration and velocity 
always are proportional to the radius r. However, we can find some radius 
r where gN( r ) < «o f° r r < r o ~ further we will call that a 'MOND radius'. 
For r < vq the gravity is in the MOND regime and the dynamics is changed. 
The evolution of the MOND radius ro in the early Universe is discussed in 
details by Sanders [4]. 

In the standard perturbation theory, if the mean density is comparable 
with the critical density of the Universe (it is true at least for large redshifts) 
the recollapse depends very strongly on the value of the overdensity. Regions 
with mean density less than the critical density will not recollapse at all and 
vice versa. Moreover, the time of the recollapse is very sensitive to the value 
of the density. In MOND it is different because there is no critical density. 
The overdensity does not play an essential role and the recollapse is similar 
for regions of different densities. 

It is quite easy to derive the linear perturbation theory in the Newto- 
nian limit (e.g. in Kolb and Turner [5]). In MOND it is much more difficult 
because there appear nonlinear terms connected with V(f, where <p is the 
gravitational potential. However, it is known that the rare highest fluc- 
tuations in the primordial density field were nearly spherically symmetric 
[6], so as long as we are concerned with the very first bound objects in the 
Universe we may assume spherical symmetry. 

Now let us consider a spherically symmetric overdensity with density 
profile 



2.1. Collapse of a pressureless fluid in MOND. 



9N = 



*M 4 

— = i^Ggr 



(4) 



g{r) = {l + 5{r))g. 



(5) 
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Of course, the density changes with time because of expansion or recollapse, 
but the mass inside some shell i 

Mi = J 4Trr 2 g(r)dr (6) 
o 

remains constant, so the deceleration for a shell of radius will be equal to 

= -/<•»> ro 

where /(rj) is the MOND gravitational force and it depends on GMi/ri and 
ao. Let us drop the subscript i. If we multiply this equation by dr/dt we 
get 



d_ 

dt 



1 (dr 

2 Vrft 



and after integrating over t we obtain 

dr 



f// v -2F(r) + C, (9) 

where F'(r) = f(r) and C is some constant which may be easily calculated if 
we know initial radius and velocity for a given function F(r). This formula 
may be easily integrated, e.g. with the Runge-Kutta algorithm. 

To show the difference between the MOND and the Newtonian gravity, 
we have performed a set of calculations. We have traced the evolution of a 
dust cloud (no pressure) expanding homogenously in all directions. Initial 
radius rj, velocity Vi, density Qi and time ti were taken from the Q = 1 
Friedmann model with no radiation for z + 1 = 500 and h = 0.65. Initial 
radius was chosen to be the MOND radius at that time. We have performed 
eight runs: 



• pure Newtonian gravity with the density equal to Qi, 0.80^ and 1.250, 

• MOND with the density as above 

• MOND with the ao parameter ten times greater than the standard 
value and mean density equal to Qi 



MOND with the ao parameter ten times lower than the standard value 
and mean density equal to Qi 
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The results are displayed at Fig. 1. Solid curves show the trajectories in 
the MOND model and the long-, middle- and short-dashed ones show the 
trajectories for the Newtonian gravity, MOND with large oo and MOND 
with low ao, respectively. For the "standard" MOND and the Newtonian 
gravity lower curves show the trajectories for runs with greater densities 
(1.25f)j) and upper curves show the trajectories for lower densities (0.80^). 
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Fig. 1. Collapse of a pressureless fluid in MOND and in the Newtonian 
gravity 

As one could expect, raising or lowering the initial density by 25% does 
not have a big influence if we apply the MOND model while it plays a crucial 
role for the Newtonian gravity if the density is near the critical one (it is 
always true in the early Universe which is flat at least asymptotically at the 
beginning) - models with higher densities tend to recollapse while models 
with the lower one do not. The results are quite sensitive to the value of 
do because it is the limit between the Newtonian and the MOND regimes - 
however, it is a quantitative effect only and the dust clouds in such models 
will always recollapse. 

2.2. Collapse of a perfect gas. 

If we take a perfect gas instead of the pressureless fluid, the evolution 
will look different because effects of pressure will moderate the recollapse, 
especially for small systems. As we assume spherical symmetry, we use 
Lagrangian coordinates. 
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The dynamics is governed by the following equations: 



dM 

dr 
dr 

dt 
dv 

dt 
du 

lit 




v 




p_dg A 

g 2 dt g' 



GM(r) 



(10) 



(12) 



(11) 



(13) 



where r is the radius of a sphere of mass M, u is the internal energy per 
unit mass, p is the pressure and g is the mass density. Here eq.(10) is the 
continuity equation, (11) and (12) give the acceleration and (13) accounts 
for the energy conservation. The last term in the eq.(13) describes cool- 
ing/heating of the gas, with A being the energy absorption (emission) rate 
per unit volume, given in details in [7]. 

We use the equation of state of the perfect gas 



where 7 = 5/3, as the primordial baryonic matter after recombination is 
assumed to be composed mainly of monoatomic hydrogen and helium, with 
the fraction of molecular hydrogen H2 always less than 10~ 3 . 

In case of modified gravity, equation (12) will look a bit different: 



where f(x) is the function inverse to the function fi(x) mentioned before 
and it is asymptotically equal to x for x clq and yja^x for x <C clq. 




In the simulations we have used the code described in [7], based on the 
codes described by Thoul and Weinberg [8] and Haiman, Thoul and Loeb 
[9]. This is a standard, one-dimensional, second-order accurate Lagrangian 
finite-difference scheme. The only changes were modification of gravity and 
putting the dark matter fraction Sldm equal to zero. However, it was neces- 
sary to make significant changes in initial conditions. 

First of all, we start our calculations at the end of the radiation-dominated 
era. For $7^ = VL m = 0.02, z eq = 203 and for fi& = 0.05, z eq = 507 as given by 
the formula provided by Hu and Eiseinstein [10], z eq = 2.50 x 10 4 r2o^ 2 ©z7i 



V = (7 - l )QU, 



(14) 




(15) 




3. Code used in the simulations 
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where 9 2 . 7 = T y /2.7K, assuming h = 0.65 and T 7 =2.7277 K. We have as- 
sumed that, as in the standard cosmology, initial overdensities may grow 
only in the matter-dominated era. We have developed and tested our own 
code to calculate the initial chemical composition and initial gas temper- 
ature. We have compared our results with the results by Galli and Palla 
[11] and they turned out to be very similar - the agreement is at a level of 
10-20%. The difference was probably due to the fact that Galli and Palla 
have included more species (e.g. deuterium and lithium) and some reaction 
rates that they used were a bit different. 



4. Results 

We have performed eight runs, for various combinations of Cl b (0.02 and 
0.05), 0\ (1 - n b and 0) and a (1.2 x 10- 8 cm/s 2 and 1.2 x 10- 9 cm/s 2 ). 
Total mass of a cloud is about 3 x 10 5 M Q . Results are displayed on Figs. 
2-7. 

Fig. 2 shows trajectories of shells enclosing 7%, 17%, 27% ... 97% of 
the total mass for the models with the "standard" value of ao- 

Fig. 3 shows the same for the models with ao 10 times lower than the 
"standard" value. 

In figures 2-5 upper plots are for = 0.02, lower ones for Qb = 0.05, 
plots at left for $7 a = 1 — and plots at right for Qa = 0. 
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Fig. 2. Shell trajectories for the first four runs, for the "standard" value of 
a . 
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Fig. 3. Shell trajectories for the last four runs, for a$ 10 times lower than 
the "standard" value. 



These results show that: 



• if we set the cosmological constant to zero, it affects the results very 
slightly - collapse is a bit slower (it is due to the fact that non-zero 
cosmological constant affects behaviour of the scale factor and, thus, 
the initial velocities) 

• time of collapse depends very strongly on ao 

• results are more sensitive to the value of for higher ao 

• mass ~ 10 5 Mq is the boundary between the regime when the gas 
effects dominate and the cloud virializes without further collapse to 
a compact object (or at least this collapse is much slower), and the 
collapse regime when the gas effects may be neglected. The boundary 
mass is a bit higher for lower ao- 

The last point needs some comments. Scientists who explore the origin 
of the Large Scale Structure often neglect gas terms, i.e. set pressure and 
internal energy equal to zero. Our paper is not the first one that deals 
with the structure formation in MOND - this problem was also explored 
by Sanders ([4] and [12]) and Nusser [13] but they were interested mainly 
in the large structure formation and they did not include gas effects. It is 
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a good approximation for large scales but gas effects play a crucial role in 
small scale structure formation. 



If we look e.g. at Fig. 2, for the outer shells the collapse is almost 
"symmetric" to the expansion, the slow-down by gas pressure is very tiny. 
In contrast, the most innermost shells show a big asymmetry - expansion, 
collapse and then shock waves make the matter to virialize (the kinetic 
energy of a collapsing shell is turned into the internal energy of the gas) so 
the collapse is stopped. Only pressure of the outer shells and the cooling 
processes (especially the cooling by molecular hydrogen H2) force them to 
collapse. 



Figs. 4 and 5 show the same as Figs. 2 and 3, but as a function of 
redshift instead of the cosmic time. They show that the collapse of the first 
compact objects in MOND is really very fast - for high fi& and "standard" 
ao it happens at z ~ 160, for lower f2& it happens at z ~ 110 and for lower 
ao at z ~ 50 — 60. It is a very strong prediction and hopefully it may be 
tested in the future by the Next Generation Space Telescope (NGST) or 
(perhaps) some its successors - if they discover luminous objects or ionized 
gas at z > 20 — 30 it will strongly support the MOND model. 



£l b = 0.02, Q A = 0.98, a =1 .2 X 1 0" B cm/s 2 Q b = 0.02, Q A = 0.00, a =1 .2 X 1 0" 8 cm/s 2 
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Fig. 4. Shell trajectories as a function of redshift, for the first four runs or 
the "standard" value of ao- 
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Fig. 6. Shell temperatures for the runs with non-zero cosmological constant 
as a function of cosmic time. 
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Fig. 7. The evolution of chemical composition for a shell enclosing 12% 
total mass, for the runs with non-zero cosmological constant. 



The last two figures show the evolution of shell temperatures (Fig. 6) 
and chemical composition of a shell enclosing 12% of the total mass (Fig. 
7) for the runs with non-zero cosmological constant. Upper plots are for the 
"standard" ao, lower ones for ao/10, plots at left for = 0.02 and plots 
at right for Qf, = 0.05. Shell temperatures behave in a very similar way: 
they fall (due to the Hubble expansion), then rapidly grow (collapse) and 
again fall (due to the cooling processes). For the "standard" ao maximal 
temperatures are significantly higher (up to 4500 K for the high-Of, model) 
so the cooling is much more efficient and the final collapse is faster. It is 
the effect of higher gravity and, thus, higher outer pressure. The evolution 
of chemical composition is qualitatively similar for all runs and all shells. It 
behaves similarly to the predictions of background cosmological models up 
to the time of collapse. The greatest differences are for and especially 
H2 which eventually reaches the final abundance of order of 10~ 3 . One 
should stress that at the time of collapse chemical reactions are much faster 
because their rates are proportional to g 2 . 

It is worth to note that in order to avoid arbitrarily small timesteps, we 
have adopted a numerical trick described by Thoul and Weinberg [8]: if a 
shell falls below some radius r c it is set to some constant value and chemical 
reactions and cooling are frozen. It means that the "flat" parts of the curves 
at the end of each plot sometimes preceeded by a strange behaviour of the 
shell temperature (very narrow "peaks" for the low-ao models) are artificial. 
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5. Conclusions 

If the assumptions of the MOND were correct, the "Dark Ages" end very 
early - about z ~ 110 — 160 for the "standard" a® and about z ~ 60 for ao 
10 times lower. It may be a good test of MOND in future because in the 
standard (A)CDM scenario the first luminous objects may appear only for 
z ~ 10 — 20. However, much more work is necessary in order to understand 
MOND and, in particular, the MOND cosmology properly. 
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